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Non-linear resonance in nearly geodesic motion 
in low-mass X-ray binaries 
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Abstract 

Wc have explored the ideas that parametric resonance affects nearly geodesic motion around a black 
hole or a neutron star, and that it may be relevant to the high frequency (twin) quasi-periodic oscillations 
occurring in some low-mass X-ray binaries. We have assumed the particles or fluid elements of an accretion 
disc to be subject to an isotropic perturbation of a hypothetical but rather general form. We find that 
the parametric resonance is indeed excited close to the radius where epicyclic frequencies of radial and 
meridional oscillations are in a 2 : 3 ratio. The location and frequencies of the highest amplitude excitation 
vary with the strength of the perturbation. These results agree with actual frequency ratios of twin kHz 
QPOs that have been reported in some black hole candidates, and they may be consistent also with 
correlation of the twin peaks in Sco X-1. 

Key words: Accretion - General relativity - X-rays: binaries - X-rays: individual (Sco X-1, J1655-40, 
J1550-564) - QPOs 



1. Introduction 

Pairs of high frequencies, known as kilohcrtz quasi- 
periodic oscillations (kHz QPOs), have been detected in 
the X-ray emission of low-mass X-ray binaries (LMXBs) 
for more than 20 neutron stars and several black holes. 
It has been suggested that a non-linear resonance within 
an accretion disc in a general-relativistic space-time met- 
ric plays a role in the excitation of the two oscillations 
(Abramowicz & Kluzniak 2001, Kluzniak & Abramowicz 
2001). 1 

It is now recognized that the frequencies of the two 
peaks in the power spectrum of X-ray variability in black 
holes are in rational ratios, lui/lu2 =m:n, with m : n = 2 : 3 
for two sources, and m : n = 3 : 5 in a third source 
(Abramowicz & Kluzniak 2001; Kluzniak & Abramowicz 
2002; Rcmillard et al. 2002, Abramowicz et al. 2002a). 
These black hole QPOs are thought to have stable fre- 
quencies, and until the discovery of the puzzling rational 
ratios they have been thought to correspond to a trapped 
g-mode or c-mode of disc oscillation in the Kerr metric 

* On leave of absence at the Center for Particle Physics, Institute 
of Physics of the Czech Academy of Sciences. 

§ On leave of absence from the Copernicus Astronomical Center 
of the Polish Academy of Sciences. 

1 Non-linear resonance has been favored as an explanation of 
QPOs also by Titarchuk (2002) in a rather different ( "transition 
layer") model. In another context, a tidally forced, Newtonian, 
non-linear vertical resonance in accretion discs has been dis- 
cussed for close binary systems (Lubow 1981; Goodman 1993). 



(Okazaki et al. 1987, Wagoner 1999, Wagoner et al. 2001, 
Kato 2001). 

In neutron star sources, the twin kHz QPO frequencies 
arc known to vary considerably and to be mutually corre- 
lated. In the context of disc oscillations, two frequencies 
varying in this fashion can be explained as the funda- 
mental and the first harmonic of the non-axisymmetric 
(m=l) g-mode (Kato 2002, 2003). In Sco X-1, a proto- 
typical Z-source, the slope of the correlation line is some- 
what steeper than 2/3, but it is not uncommon to find the 
source in a state when lo\ « 600 Hz, and cl>2 ~ 900 Hz, i.e., 
at about 2 : 3 ratio (Abramowicz et al. 2002b). 

It would appear that the 2 : 3 frequency ratio is com- 
mon to neutron star and black hole systems. Kluzniak & 
Abramowicz (2002), and Abramowicz & Kluzniak (2003) 
suggest that the 2 : 3 ratio follows from the properties 
of parametric resonance between two eigenfrequencies, of 
which one is always lower than the other. Specifically, 
two modes of disc oscillation arc thought possible with 
frequencies close to the two epicyclic frequencies of free 
orbital motion. Within general relativity, the radial 
epicyclic frequency is generically lower than the merid- 
ional one, their ratio varying from one to zero as the ra- 
dius of circular orbits decreases from infinity down to that 
of the marginally stable orbit. 

The observed rational ratios of frequencies may follow 
from fundamental properties of both strong gravity and 
a non-linear resonance between radial and vertical oscil- 
lations in accretion discs with fluid lines that are nearly 
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Fig. 1. Instability regions of the Mathicu equation (8) 
with the amplitude of variation of the eigenfrequency hi, 
n = 2u)g/uj r , and A = 0. Levels of shading indicate different 
values of the growth rate of 58(t). Three tongues of insta- 
bility, of the order n = 1, 2, and 3, are clearly visible. The 
equations (2)-(3) considered in this paper have similar regions 
of instability. Note that the range of frequencies for which the 
instability develops increases with hi. 

geodesic, nearly circular, and nearly planar. In this pa- 
per we investigate a simple mathematical model of such 
motions. To this aim, we discuss the properties of para- 
metric resonance in the Paczyhski-Wiita (1980) model of 
the Schwarzschild metric. 

2. Nearly-geodesic motion 

Let us consider a fluid in which the flow lines are only 
slightly non-circular, and slightly off the # = 7r/2 symmetry 
plane. In spherical coordinates, 

r(t)=r + Sr(t), 6(t) = ±tt + 60(t), (f>{t) = Slt. (1) 

With accuracy to the third order in 5r <C an d 56 <C n/2, 
equations of fluid motion read: 
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Here the time derivative is denoted by a dot, the ra- 
dial derivative by a prime, and fi are components of a 
small force of non-gravitational origin (pressure, viscous, 
magnetic, or other). The epicyclic eigenfrcquencics tog 
and oj r arc defined, in terms of the effective potential 



U = 3>(r,9) + £ 2 /(2r 2 sin 2 6) and the specific angular mo- 
mentum I = <j)r 2 sin 2 6, by 
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Derivation of these equations assumes equatorial plane 
symmetry for the gravitational potential <f>. Here we 
adopt a spherically symmetric potential 2 
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From eqs. (5) and (6) it follows 

2 _ GM f r ms \ 2 _ GM 

1 < LUn — 



r(r-r G ) 2 ' 



(6) 



(7) 



where r ms = 3r G is the marginally stable orbit. 

3. Strong gravity's 2 : 3 resonance 

A parametric resonance instability occurs near uj r = 
2ujg/n for n— 1,2,3,..., in an oscillator that obeys the 
Mathieu type equation of motion (Landau & Lifshitz 
1976), 



- u 2 [1 + h x cos(uvt)] 56 + Xc 



0. 



(8) 



When the coupling is weak (hi small), the strongest in- 
stability occurs for the lowest value of n possible (see 
Figure 1). 

With this textbook mathematics in mind, let us con- 
sider a simplified version of eqs. (2)-(4), with fi = 0, and 
some of the higher order terms neglected. With accuracy 
up to linear terms in (3) one obtains, obviously, a solution 
8r(t) oc cosLOrt. Substituting this solution in (2) brings this 
equation close to the standard form (8). We have solved 
these equations for nuj r = 2ujg, with n a positive real pa- 
rameter. The instability regions for a lower-order version 
of equations (2)-(3) resemble the well-known tongues of 
instability of the Mathicu equation. Since u r < ujg in gen- 
eral relativity, the lowest value for which the resonance 
can occur is n — 3, and the first two tongues of insta- 
bility in Fig. 1 will be absent; i.e., the ratio of the two 
eigenfrcquencies in (8) is 2 : 3 for the strongest resonance 
(Kluzniak & Abramowicz 2002). 

However, the behaviour of the resonance can only be 
properly studied if the third-order terms in (2) and the 
second-order terms in (3) are retained. These terms pro- 
vide the non-linearity necessary to saturate the ampli- 
tude at a finite value and the damping which affects the 
frequency of oscillations. Using the standard analytic 
method of successive approximations (Landau & Lifshitz 
1976), we verified that, as expected, no parametric reso- 
nance occurs for strictly geodesic motion, i.e., for fi = in 
eqs. (2)-(4). We also checked numerically that there is no 
resonance in the geodesic case. However, the resonance 

2 See Paczyhski & Wiita (1980). This form is known to be a 
convenient model for the external gravitational field of a non- 
rotating black hole or a neutron star; it captures the essential 
feature of Einstein's strong gravity relevant here, i.e., the fact 
that oj r (r) < utg(r), and w r = at the marginally stable orbit. 
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Fig. 2. Left: Location of the first three resonances is shown for a particular choice of initial amplitude, 89(0) =0.5, by encoding the 
growth rate of \S8\ from eq. (2) with different levels of shading in the (r, a)-plane (radius versus strength of perturbation). Note that 
a decreases along the y-axis. For a = 1.00, the resonant orbits are located at radii r m :n (indicated by vertical lines) where u) r and ujg 
are in rational ratios, m: n. In particular, the trace of the 2 : 3 resonance (n = 3) starts from ro = 9.2GM/c 2 in the Paczyhski-Wiita 
potential. The frequencies and radius of resonant orbits vary when a is changed. Right: The maximum amplitude A of S6(t) for a 
fixed value of a = 0.98 as a function of the radius of the unperturbed circular orbit. 



does occur (as anticipated in Kluzniak & Abramowicz 
2002) for even very slightly non-geodesic motion, when 
the higher order terms are influenced by non-geodesic ef- 
fects of a rather general form and right amplitude. 

4. Non-geodesic coupling 

In general it is not possible to specify the form of the 
very small force fi exactly, because the nature of non- 
geodesic forces in accretion discs is not yet very accurately 
known. Here, we only explore the basic mathematical 
form of the solutions, not the physical origin of the fi 
force. We assume ad hoc that the possible non-geodesic 
effects have the form of a non-linear isotropic coupling 
between the two components of geodesic deviation: 



fe = +a e 



fr 
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(9) 



(10) 



Here, a r and ag measure the strength of non-geodesic 
forces in comparison with the geodesic terms of corre- 
sponding form and order. We first assumed a r = ag = a 
for the sake of simplicity. 3 

After choosing a particular value of a, we solved eqs. 
(2)-(3) numerically by using the Runge-Kutta adaptive 

3 Perturbations considered here may be consistent with the for- 
mation of non-axially symmetric coherent structures in accre- 
tion discs ( "planets" , "vortices" , "magnetic flux tubes" or "spi- 
ral waves"), as proposed by several authors: Abramowicz et 
al. (1992); Goodman, Narayan, & Goldreich (1987); Adams & 
Watkins (1995); Bracco ct al. (1998); Karas (1999); Kato (2002, 
2003). 



step-size routine. The integration was performed over a 
wide range of starting radii, 6 < ro(c 2 /GM) < 15. For 
the initial conditions, we took 56(0) = 0.5, 59(0) = 0.01, 
5r(0) = 0.5, and 5r(0) = 0.01. We then repeated the same 
procedure for another value of a, and in this way we con- 
structed a sequence of solutions that is parameterized by 
a and ro. 

One can locate radii r p (a) for which \59(t;r p ;a)\ grows 
in time up to a maximum amplitude. The saturation am- 
plitude, A, is shown in Fig. 2 (right panel) func- 
tion of radius. Clearly, the largest amplitudes correspond 
to parametric resonance. Close to those radii the test- 
particle epicyclic frequencies are in rational ratios n : m; 
r P = r 3:2: r 4:2i an d r^. Note that the growth rate of 
the amplitudes |<W(i)|, \5r(t)\ decreases with decreasing 
a. The 2 : 3 resonance is the strongest one, as expected. 

The exact slope of the 2 : 3 resonance trace r p (a) , start- 
ing from the point (ro,a) = (9.2,1), depends on the choice 
of coupling parameters in eqs. (9)-(10), so the above as- 
sumption of identical values for ag = a r is not univer- 
sal. For example, assuming a r = nag and setting k = 0.2, 
one finds that the end of the trace is moved to the point 
(9.0,0.95). It is thus shifted in the (r, a)-plane with re- 
spect to the case n = 1, for which the trace ends at the 
point (8.8,0.95), as shown in the left panel of Fig. 2. Also 
the initial conditions for 58 and 5r have an influence on 
the trace slope. Therefore, the corresponding frequencies 
that are present in the solutions for 5r(t), 59(t) depend on 
the details of the solution. Nevertheless, the ratio of the 
frequencies stays at 2 : 3 at the maximum amplitude. This 
fact is in agreement with the intuition about the proper- 
ties of eqs. (2)-(3), and we checked it by Fourier-analyzing 
their solution. 
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Fig. 3. Pairs of oscillation frequencies at points along the 
n = 3 trace of Fig. 2 were obtained by integrating eqs. (2)— (3) 
with a different initial value 5(0) = 0.8, yielding a frequency 
ratio differing from 2:3. The solid line is a result of the 
calculation described in the text. For comparison with ob- 
servations, we have scaled all computed frequencies with one 
arbitrary multiplicative factor. The frequencies then depend 
only on the strength of perturbation a, running in the figure 
from a = 1.0 at uj* = 900 Hz to a = 0.95 at about 1050 Hz. 
Note the agreement between the slope of the computed line 
and the observed kHz QPO frequencies in Sco X-l, also shown 
(as circles). The dot-dashed line is the least-squares best-fit 
to the data points (two lines with a 3% offset from the best-fit 
line are also plotted). The dashed line gives a reference slope 
of 2/3 for comparison. The Sco X-l data (van der Klis et al. 
1997) has been kindly provided by Michiel van der Klis. 

In this paper we assumed the gravitational potential 
of eq. (6) for the sake of simplicity. When a more ac- 
curate, full formulation in general relativity is developed 
together with a physical mechanism responsible for the 
non-geodesic coupling term (9)-(10), one will be able to 
constrain the absolute value of the mass of the accreting 
body. 

5. Frequency correlation in Sco X-l 

As noted in Section 1, the frequency ratio of the QPOs 
in some black hole sources is accurately 2:3. Specifically, 
in GRO J1655-40 the frequency ratio is that of 300 Hz and 
450 Hz, while in XTE J1550-564 it is that of 184 Hz and 
276 Hz (Abramowicz & Kluzniak 2001, Rcmillard et al. 
2002). The properties of parametric resonance, discussed 
above, provide a natural explanation for this ratio. 

Can the same mechanism be responsible for the ob- 
served frequencies of the kHz QPOs in neutron stars? In 
neutron star sources the two kHz QPOs vary considerably 
in frequency, and their ratio is not always 2 : 3. A case in 
point is Sco X-l, whose two kHz frequencies are linearly 
correlated (Fig. 3), but not directly proportional to each 
other (van der Klis et al. 1997). However, the distribu- 
tion of points along the line of correlation is not uniform, 
and the resulting distribution of frequency ratios has a 



prominent peak at about 2 : 3 value, strongly suggesting 
the presence of a non-linear resonance (Abramowicz et al. 
2002b). 

We note that QPOs are not coherent oscillations, and 
that the typical integration time (minutes) greatly exceeds 
the coherence time of the signal (fraction of a second). 
We find that the slope of Sco X-l frequency correlation 
can be reproduced if we assume that the signal observed 
at a given time corresponds to the frequencies originat- 
ing at a particular radius, for example the radius r p (a) 
discussed in the previous section. If this is not strictly 
equal to the radius where the maximum amplitude is 
attained, the ratio of the frequencies departs from 2/3. 
Specifically, taking an arbitrary initial amplitude 69(0) 
we constructed Fourier power spectra of S9(t;r r , a), and 
Sr (t ; r p , a) . These power spectra peak at frequencies ujg (a) 
and to* (a), which we interpret here as the observed pair 
of oscillations. 

Figure 3 shows the correlation between ujg and to* along 
the a-trace from Fig. 2. We assumed k = 1 as before. 
The actual value of 58(0) has been selected to match the 
slope of Sco X-l data (van der Klis 1997). Because a 
measures the deviation from geodesicity, its value is likely 
to be a function of the accretion rate in the disc, a(M). 
This notion agrees with the reported correlation between 
the frequencies and the position of the source along the 
Z curve in the color-color diagram. Fig. 3 illustrates that 
the simple scheme that we introduced in the present paper 
may, indeed, explain the slope of u>*(ojg) relation. 

6. Conclusions 

We have expanded, through the third order, relativistic 
equations of test motion about a circular geodesic, and 
have shown that the lowest order expansion corresponds 
to Mathieu's equation with non-standard damping. As 
expected, the solutions of these (unperturbed) equations 
do not show any resonant phenomenon — they simply de- 
scribe geodesic motion. 

When a sufficiently large non-geodesic perturbation is 
imposed, the deviations grow rapidly at certain resonant 
radii and the motion departs from circular geodesies. This 
is caused by parametric resonance between the meridional 
and radial epicyclic motions, and the strongest resonance 
occurs when the two dominant frequencies are near ratio 
2 : 3. 
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